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ABSTRACT 


Homotopy continuation, a powerful tool for identifying steady state multiplicities, 
is applied to reactive distillation systems to assess the impact of these multiplicities on 
column design, operability and control. The technique is applied in two distinct ways. In 
the first, continuation with respect to the catalyst weight is used to analyze multiplicities 
as reaction is introduced in simple distillation. An MTBE and a Methyl Acetate RD 
column are used as examples. For the MTBE column, upto five steady states in the 
kinetically controlled regime are detected. Of these, only three remain as the catalyst 
weight is increased and reaction equilibrium is approached on the trays. For the methyl 
acetate column, three steady states in both the kinetically controlled and reaction 
equilibrium regimes are obtained at fixed reflux rate and reboiler duty. If the column 
specification is changed to fixed reflux ratio and reboiler duty, only a single steady state 
solution is obtained for all catalyst weights. The result shows the substantial impact of 
column operating policy on the existence of multiple steady states. 

In the second application of homotopy continuation to RD columns, continuation 
with respect to an input column specification is used to study its effect on output 
variables such as tray temperatures and compositions. The input-output relationships are 
systematically analyzed to synthesize robust column operating structures that ensure high 
product purity and reaction conversion in the face of major disturbances. The basic 
philosophy is to prefer control structures with nearly linear input-output relationships and 
avoid steady state multiplicities. Sensitivity analysis is used to obtain the sensitive tray 
location for the temperature / composition sensors used in a control structure. An MTBE 
RD column is used as an illustrative example. 
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Chapter 1 


Introduction 

Reactive distillation (RD) combines reaction and separation in a single column. It 
is a promising technology for process intensification with significant economic benefits 
over conventional “reactor followed by separator” processes. Besides the rectifying and 
stripping sections, reactive distillation columns comprise a reactive section in which 
reactants are converted into products. A typical RD column is shown in Figure 1.1. 
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Figure 1.1: A t)^ical reactive distillation column 
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Commercial application of this process is mainly in esterification systems such as 
methyl acetate and etherification systems such as methyl tertiary butyl ether (MTBE), 
ethyl tertiary butyl ether (ETBE), and tertiary amyl methyl ether (TAME). Methyl acetate 
is a general solvent used for cellulose, lacquers, paints, resins, coatings, and perfumes and 
is also used as an intermediate in the manufacture of pharmaceuticals, artificial leather 
and synthetic flavoring. MTBE is valuable as gasoline additive for increasing the octane 
rating of the fuel. It is known as an "oxygenate" because it raises the oxygen content of 
gasoline and helps in its complete combustion reducing harmful CO emissions. MTBE is 
the key octane enhancing ingredient of unleaded gasoline. It has replaced the eco- 
unfriendly tetra ethyl lead (TEL), the octane enhancer in leaded fuels. ETBE and TAME 
are alternatives to MTBE for enhancing the fuel octane rating. 

Reactive distillation when feasible has several advantages over conventional 
reactor-separator process: 

1. Nearly complete conversion in equilibrium limited reactions due to the continuous 
removal of reaction products from the reactive zone. This is a consequence of the Le 
Chetalier’s principle. Complete (increased) conversion eliminates (reduces) the costs 
associated with recycle. 

2. Enhanced selectivity of the main reaction due to the removal of one (or more) 
reactants in the side reactions by volatilization. 

3. Better heat utilization in the case of exothermic reactions where the heat of reaction is 
directly utilized for fractionation, reducing the reboiler duty. 

4. Elimination of azeotropes due to reaction with other components .The azeotrope is 
literally “reacted away”. 
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5. Process intensification with reaction and separation taking place in a single unit 
operation. 

The methyl acetate reactive distillation process, first commercialized by the 
Eastman Chemical Company, is a classic success story in reactive distillation history. The 
improvement over the commercial process is so dramatic that the capital and operating 
costs reduce to one fifth of the conventional process. Figure 1 .2 shows a schematic of the 
two processes. The main is as follows: 

CH.COOH + CH,OH ^=± CHfOOCH, +H^O 
(HOAc) ( MeOH) {MeOAc) 

The traditional process consists of one reactor and a train of nine distillation 
columns. The separation of the reactor effluent into pure methyl acetate is complex due to 
the presence of azeotropes between methyl acetate and methanol, methyl acetate and 
water, and a near azeotrope between acetic acid and water. The RD process, in contrast, 
is a single column producing methyl acetate of greater than 95% purity with nearly 100% 
conversion. The column is thus a complete plant in itself 

Another commercially important example is MTBE production using RD 
technology. The production of MTBE is also traditionally carried out in a reactor 
followed by distillation column. MTBE is derived by the reaction of isobutene and 
methanol. The conventional and RD processes are shown in Figure 1.3. The main 
reaction in MTBE production is: 

{CH^\C = CH^ + CH^OH ^=± {CH,),COCH^ 

{Isobutene) {MeOH) {MTBE) 
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(a) Conventional process (b) Reactive Distillation 


Figure 1 .2: A pictorial comparision between the conventional process and reactive 
distillation for methyl acetate production. Adapted from Siirola (1995) 
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Figure. 1 .3: MTBE production by conventional and reactive distillation processes 



It is emphasized that for RD to be a feasible production scheme, the relative 
volatility of the components should favor high reactant concentration in the reaction 
zone. Also the reaction rates at the mixture bubble point temperature should be large. RD 
is thus not an alternative for the production of each and every chemical. However, 
whenever feasible, the economic advantages of RD over conventional processes can be 
phenomenal. 

The economic advantages of RD not withstanding, the operation and control of 
RD columns is a complicated task that can be a hindrance in the successful 
implementation of the technology. The coupling of reaction and separation in a single 
unit makes the system highly nonlinear so that predicting column behavior in the 
presence of disturbances and changed operating conditions is non-trivial. Consequently, 
formulating the best operating strategy can be quite complex and non-intuitive and the 
strategy may change depending on the operating conditions. The issue is further 
complicated by the existence of multiple steady states that are to quite common in RD 
systems. A complete and thorough understanding of these multiplicities is a must for 
proposing robust control strategies for the proper regulation of the column in the face of 
disturbances. 

We distinguish between two types of steady state multiplicities: input and output 
multiplicities. The former refers to different inputs (standard column specifications such 
as reboiier duty or reflux ratio) giving the same output (conversion or purity or distillate 
rate) while the latter refers to the same input giving different outputs. The idea of input 
and output multiplicity is illustrated in Figure 1.4(a) and 1.4(b) respectively. 
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OUTPUT VARIABLE 


In the literature on the simulation of RD columns, conventional algorithms such 
as Naphtali - Sandholm and Inside - Outside methods have been adopted for solving the 
governing MESH equations. These methods, depending on the initial guess can converge 
to any one of the multiple steady states, in case multiple steady states do exist. The 
methods are thus not geared towards detecting steady state multiplicities. A systematic 
approach for studying the steady state multiplicities is needed and this forms the main 
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Figure 1 .4: Input and output multiplicities 

focus of this work. As shown in Figure 1.4, the steady state multiplicities are clearly 
revealed in the solution diagram. Traditional equation solvers give a point on the solution 
diagram (I, II or III) depending on the initial solution guess. Given a solution, the input 
(or output) can be varied about this solution to trace the complete solution path. This is 
refeiTed to continuation with respect to a parameter. Thus continuation methods are a 
natural choice for obtaining the solution diagram and revealing multiplicities. 
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In order to start the continuation, one point on the solution branch is needed. This 
initial solution can be obtained in two ways as: a) the steady state solution to an ordinary 
distillation problem with no reaction and then introducing reaction as the continuation 
parameter and b) the steady state solution of the RD problem from conventional 
distillation solvers such as the Naphtali-Sandholm or inside-outside method and then 
varying an input (or output) about that steady state. In this work, RD columns of 
commercial interest viz Methyl Acetate and MTBE are studied for steady state 
multiplicities using approach (a). Approach (b) is applied for formulating the control 
strategy for a MTBE RD column. 

Continuation using approach (a) described previously, is implemented using the 
catalyst weight as the continuation parameter. The ordinary distillation (no reaction) 
solution is obtained at zero catalyst weight and the desired reactive solution is 
approached by slowly increasing the catalyst weight. The homotopy path is thus traced 
from zero reaction extent to the design reaction extent. In the MTBE system, a reasonably 
large zone with five multiple steady states is obtained vath the high conversion steady 
state corresponding to the design steady state. As the catalyst weight is further increased, 
two of these steady states disappear. Thus there are five steady states in the kinetically 
controlled RD regime in contrast to only three steady states in the reaction equilibrium 
controlled RD regime. For the methyl acetate column, three steady states in both the 
kinetically controlled and reaction equilibrium regimes are obtained at fixed reflux rate 
and reboiler duty. If specification set was changed to fixed reflux ratio and reboiler duty 
then no steady state multiplicities are seen. .Using approach (b) as noted above with the 
high conversion steady state obtained using the NS method as the initial solution, input 
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and output multiplicities are clearly seen in MTBE system. Depending on the input 
variable (reboiler duty, feed flow rates) used for the continuation, there can be multiple 
solutions for the output variables. These results have important implications on 
formulating the control strategy for a RD column. Input variables that lead to nearly 
linear relationships with the output variables are preferred as potential manipulated 
variables. When output multiplicities cannot be avoided, inputs that have the design 
steady state far away from the region of multiplicities are preferred. The power of this 
simple rationale in the synthesis of RD column control structure screening is illustrated 
for the MTBE system. 

The thesis is organized as follows. In Chapter 2, the literature related to multiple 
steady states in reactive distillation with emphasis on MTBE and methyl acetate is 
surveyed. A brief overview of the literature on homotopy continuation methods is also 
presented. Chapter 3 describes the homotopy continuation method in detail. In particular, 
the concept of dynamic reparametrization for tracing the solution diagram around turning 
points is discussed. The homotopy continuation method using approach (a) is applied to 
the methyl acetate and MTBE systems in Chapter 4. Homotopy continuation using 
approach (b) is applied to MTBE system for control structure synthesis in Chapter 5. The 
implications of the observed input / output multiplicities on control structure synthesis are 
discussed. Finally Chapter 6 presents the conclusions that can be drawn from this work 
and provides recommendations for future work. 
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Chapter 2 


Literature Review 


This chapter briefly reviews the literature on reactive distillation and homotopy 
continuation. These distinct topics are reviewed in separate sections. Section 2.1 focuses 
on reactive distillation literature. In particular, emphasis is laid on reviewing published 
work on input output multiplicities and multiple steady states in reactive distillation 
systems. Section 2.2 provides a brief overview of the literature on homotopy continuation 
method. In particular, applications of the technique to RD columns are highlighted. 

2.1 Reactive Distillation Literature 

Reactive Distillation (RD) is an old process that combines reaction and separation 
in a single column. The first patent dates back to the early 1920s and was awarded to 
Backhaus for esterification systems (Backhaus (1921), (1922) and (1923 a, b)). The 
Eastman Kodak Company developed the first commercial RD process for the production 
of methyl acetate via the esterification of acetic acid and methanol with H 2 SO 4 as the 
catalyst (Agreda and Partin (1984)). This process is a classic success of RD technology 
with reported capital and operating costs being a fifth of the traditional process (Agreda, 
(1990)). Given its potential, extensive research and development work in the past two 
decades have lead to commercial processes for MTBE, ETBE and TAME production 
(Degarmo (1992)). 
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RD can be advantageous in driving the conversion of equilibrium limited 
reactions to near completion due to the removal of products by vaporization, heat 
integration with the reaction heat being used for volatilization, side reaction suppression 
by the removal of side reaction precursors from the reaction phase and “reacting away” of 
azeotropes (Malone and Doherty (2000)). The realization of all these advantages requires 
favorable vapor liquid equilibrium and substantial reaction rates at the colunm operating 
temperatures as dictated by the tray mixture bubble point temperature. When realized, 
these advantages lead to process intensification with the complete plant being 
miniaturized into a single RD column, as in the case of methyl acetate production. 
Significant economic savings are therefore realized. Luyben (2004) has recently 
performed a quantitative comparative study on the economic benefits of RD over a 
conventional “reactor followed by separator” process. Ideal VLB is assumed in the study. 
The paper shows that when the relative volatilities of the components are favorable, the 
RD process is cheaper than the conventional process by a factor of 2-3 while the 
conventional process is the preferred option in case the relative volatilities are 
unfavorable. 

The economic advantages provide a compelling reason for adopting RD as a 
technology. However, the operation and control of RD columns a particularly challenging 
task. Compared to conventional processes, the control degree of freedom is reduced in 
RD processes so that the separation and reaction extent must be regulated using fewer 
available valves. Also, the coupling of reaction and separation leads to high non-linearity 
in the governing equations leading to multiple steady states and non-monotonic 
temperature profiles. The reader is referred to Taylor and Krishna (2000) for a 
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comprehensive review on the modeling of RD columns. A thorough understanding of RD 
column behavior through extensive modeling and simulation is a must for devising robust 
control strategies for the stable and economic operation of the RD columns in the face of 
disturbances. For this purpose systematically conducted simulations provide vital insight 
on the governing input-output relationships. 

In this context, several researchers have, over the years, come across multiple 
steady states in RD column simulations. There are several reports oii multiple steady 
states in RD literature. After the initial skepticism, multiple steady states have come to be 
accepted as routine phenomena in RD systems. Assuming the equilibrium stage model, 
Jacobs and Krishna (1993) reported a low conversion and a high conversion steady state 
in an MTBE column. Nijhuis et al (1993) have also reported multiple steady states in an 
MTBE column. Huan et al (1994) attempted to explain these multiplicities as due to the 
substantial increase in the activity of methanol in the presence of the inert n-butane 
leading to lower temperatures in the column with consequent lower reaction conversion. 
Chen et al (2002) used arc-length continuation to obtain upto three steady states for the 
MTBE system. They also reported three steady states for the TAME system in the 
kinetically controlled regime with only one steady state in the limit of reaction 
equilibrium on all the reactive stages. Higler et al (1999) reported multiple steady states 
in MTBE RD using the non-equilibrium stage model. Guttinger and Morari (1999a and 
1999b) performed an oo/oo analysis, where the two infinities refer to total reflux. and 
infinite number of stages, for predicting multiple steady states in RD columns. They use 
the MTBE system as an example case study in their work. Ciric and Miao (1994) 
explored steady state multiplicities in an RD column producing ethylene glycol from 
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ethylene oxide and water. They observed two sets of multiplicities, a three branch 
multiplicity that occurs at large reactive stage holdup and a complex switchbacking 
multiplicity that occurs at small reactive stage hold up. 

All the above articles on multiple steady states are based on simulations only. 
Reports detailing much needed experimental evidence have also been published in the 
open literature. Rapmund et al (1998) have experimentally validated the existence of 
multiple steady states in the production of TAME. Mohl et al (1999) have experimentally 
validated the existence of multiple steady states in both the MTBE and TAME systems. 

As mentioned earlier, the existence of multiple steady states and input-output 
multiplicities make the task of RD column regulation using appropriate control strategies, 
a challenging one. In recent years, many articles have appeared in the literature on the 
operation and control of RD columns. Prominent among these is the pioneering work of 
Luyben (1999, 2000, 2002, 2002 (a) and 2002 (b)) who have published a series of articles 
on the control of RD columns. Their work offers much insight on the key issues that 
govern the synthesis of control structures for RD columns. Control of an RD column with 
ideal VLE is studied in Al Arfaj and Luyben (2002). The authors show that controlling a 
tray temperature in the stripping section and a composition in the reactive section 
provides good control for production rate changes. The control of the industrially 
important methyl acetate column with highly non-linear VLE is studied in Al Arfaj and 
Luyben (1990). A control structure that maintains a tray temperature in the reactive 
section and the stripping section using the acetic acid feed and methanol feed 
respectively, with fixed reflux ratio and the reboiler heat input as the production rate 
handle provides the best performance of all the control structures studied. It is worth 
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mentioning that this is the control system implemented by Eastman for operating their 
column. The results show that composition measurement based control structures that 
work well for the ideal column do not work for the methyl acetate column. This is 
attributed to the high nonlinearity between the input output relationships for the 
composition based control structures. The non-linearity is much less severe if temperature 
is used as an inferential measurement so that good control is obtained. A1 Arfaj and 
Luyben have also studied the control of an olefm-metathesis column (A1 Arfaj and 
Luyben (1999)), an ethylene glycol column (A1 Arfaj and Luyben (2000)) and the design 
and control of a RD based TAME process (A1 Arfaj and Luyben (2002)). 

Other authors have also studied the operation and control of RD columns. These 
include N Vora and Doutidis (2000) who study the non-linear control of an ethyl acetate 
column. Sneesby and ‘O Tade (1999) study the control of an ethylene glycol column. 
Pattern based control of an ETBE column is studied by Y.-C. Tian et al (2003). Most of 
the work is however esoteric and un-necessarily mathematical and does not provide much 
intuitive insight on the operation and control of RD columns. A notable exception is the 
work of Wang, Wong and Lee (2003) on the control of an MTBE column. They show 
that good control is obtained when the reflux ratio and the feed ratio is kept fixed and a 
tray temperature in the stripping section is maintained using reboiler heat duty. A tray 
composition must be controlled to maintain the stoichiometric feed balance by 
manipulating the isobutylene feed in a cascade arrangement. The fresh methanol feed 
thus acts as the production rate handle. 

The key -to devising good control structures for RD columns is to choose the 
manipulated (input) and controlled variables (output) carefully so that the input-output 
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pairings avoid input or output multiplicities. This idea is well illustrated in the paper by 
Wong, Wang and Lee. Given that these input-output relationships can be quite complex 
for RD columns with turning points causing multiplicities (see Figure 1.4), homotopy 
continuation about the design operating conditions is used to obtain these relationships. 
This is done to avoid convergence problems near the turning points. A brief review of 
homotopy continuation as applied to ChE systems (RD in particular) is provided in the 
next section. It is noted that in this work, we also use homotopy continuation with respect 
to catalyst hold up in a reactive stage to analyze for multiple steady states in RD systems. 

2.2 Homotopy Continuation Literature 

Steady state modeling problems in chemical engineering involve the solution of 
nonlinear algebraic equations. This is true for isolated unit operations such as reactors 
and distillation columns as well as for complete plants with material and energy recycle. 
In case of highly non-linear equations, convergence using traditional methods such as the 
Newton Raphson method may be difficult. More importantly, in case of multiple steady 
states, the locally convergent Newton Raphson method would converge to the nearest 
solution depending on the initial guess provided and miss the other solutions. In such 
situations, homotopy continuation methods are a possible alternative for obtaining all the 
steady state solutions to the governing equations. 

The basic idea behind homotopy continuation is to start from a known steady state 
solution and then slowly vary a parameter about this steady state and trace the complete 
solution diagram (Pushpavanam (2000); and Kubicek and Marek (1977)). In cases where 
multiple solutions exist, turning points that are bifurcations of the simplest type, may be 
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encountered so that the solution branch turns back around the turning point. This is 
illustrated in Figure 1 .4. The existence of multiple steady states for a given value of the 
continuation parameter is thus clearly revealed in the solution diagram. 

Homotopy continuation may also be used to solve a difficult problem g(x) = 0 by 
first starting with an easy problem with known solution f(x) = 0. A mixed homotopy 
function h(x) is then defined by mixing f(x) and g(x) using a mixing parameter t. The 
homotopy function is then solved by continuing along the mixing parameter t from the 
limit where h(x) = f(x) for which the solution is known, to the limit where h(x) = g(x). In 
some sense, the simple function f(x) is bent or twisted slowly to the complex function 
g(x) using continuation along the parameter t. An example for the homotopy function is 
h(x) = t.g(x) + (l-t).f(x). Thus t is varied from 0 to 1 with h(x) at t = 1 corresponding to 
the desired solution of g(x) = 0. 

JD Seader at the University of Utah and WD Seider at Penn State have pioneered 
the application of homotopy continuation to chemical engineering problems. A good 
review of the method is available in Waybum and Seader (1987) and the interested reader 
is referred to the same. Lin et al apply homotopy continuation for computing multiple 
solutions to systems of interlinked separation columns. Kovach and Seider, (1987) apply 
homotopy continuation methods to azeotropic distillation columns. The algorithm of 
Allgower and Georg (1981) is applied to obtain two steady states for an industrial 
azeotropic distillation column over a very narrow range of reflux ratios. In addition; the 
authors detect five steady states for the classic ethanol-water-benzene system for the 
dehydration of ethanol. Jalali and Seader (1999) have used homotopy continuation for 
robust multi-phase and multi-reaction equilibrium calculations. Another interesting 
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of method is for calcuiating all the azeotropes in a general 
multicomponent vapor-liquid mixture (Malone and Doherty (1998)). The algorithm is 
used in the commercial process simulator Aspen Plus. 

Homotopy continuation methods have been extensively used in the literature for 
the simulation of reactive distillation systems. Chang and Seader (1987) use homotopy 
continuation for solving the highly non-linear governing MESH equations for the steady 
state simulation of an ethyl acetate RD column for different sets of specifications. 
Sneesby et al (1998) apply homotopy continuation to obtain multiplicities in ETBE and 
MTBE RD columns. A high conversion and a low conversion steady state are found in 
the solution diagram for the ETBE system. In the MTBE system, the authors have 
reported a narrow operating zone with five steady states at constant reboiler duty. Chen et 
al perform a bifurcation study using homotopy continuation with the Damkohler number 
as the continuation parameter, to obtain zones with three steady states in TAME and 
MTBE systems. Ciric and Miao (1994), apply homotopy continuation for the detection of 
multiple steady states in an ethylene glycol RD column using catalyst weight as the 
continuation parameter. 

In addition to its widespread use for the detection of multiple steady states in RD 
systems, homotopy continuation has also- been used for understanding the complex input 
output relations in RD systems. Al Arfaj and Luyben (2002) use homotopy continuation 
to detect input and output multiplicities in a methyl acetate RD column. Plots between the 
distillate rate and the reboiler duty at constant reflux rate or reflux ratio show the 
operating zones of output multiplicity. These plots clearly show that operating at constant 
reflux ratio is preferable over constant reflux as the zone of output multiplicity is very 


17 



narrow in case of the former. In a similar vein, Wang et al (2003) study steady state input 
output relations for an MTBE column to propose effective control strategies. 

This completes a brief review of the literature on reactive distillation and 
homotopy continuation. Homotopy continuation will be extensively used in the our work 
for control structure screening and a systematic study of multiple steady states in 
commercially important RD systems for methyl acetate and MTBE production. 
Accordingly, the next chapter describes the mecheinics of homotopy continuation as 
applied in our work to RD systems. 
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Chapter 3 


Homotopy Continuation Methods 


The primary motivation behind this thesis is the systematic steady state analysis 
of RD columns for obtaining simple and robust control strategies. In this context, 
homotopy continuation is extensively used as a tool for understanding the steady state 
behaviour of RD columns. Continuation is applied to RD colixmns in two distinct ways. 
In the first case, continuation with respect to catalyst weight (or reactive tray hold-up) is 
used to study the presence of multiple steady states. The catalyst weight per tray is an 
important design variable that must be chosen so that steady state multiplicities are 
minimized. In the second case, continuation about a design steady state is used for 
studying the relationships between potential manipulated (input) and controlled (output) 
variables to verify the presence (or absence) of input or output multiplicities. A 
systematic study would reveal potential input output pairings that avoid multiplicities so 
that good control structures that ensure robust regulation of the column in the presence of 
disturbances are obtained. 

Given the central role of homotopy continuation as a steady state analysis tool, 
this chapter describes the theory of popular homotopy continuation algorithms and its 
formulation for RD columns as applied in this work. The principal challenge in using 
homotopy continuation for obtaining the complete solution diagram is the strategy for 
continuation around turning points. Two popular methods in the literature are dynamic 
reparameterization and arc-length continuation. These are discussed in Section j.l. 
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Section 3.2 formulates the homotopy continuation problem for RD columns for the two 
cases as discussed in the previous paragraph: (a) continuation with respect to the catalyst 
weight and (b) continuation with respect to operating conditions about a design steady 
state. 

3.1 Homotopy Continuation Algorithms 

Let us say we are interested in finding the variation in the steady state solution of 
a system of equations 

h(xA) = 0 (3.1) 

with respect to the parameter X. Assuming that xq is a known solution at A, = )to so that 
h(xo, X(j) = 0, the parameter can be changed incrementally to A, + A.X, and the 
corresponding solution xo + Ax obtained using iterative equation solving algorithms such 
as the Newton Raphson method. The iterative algorithm will converge easily as the new 
solution Xo + Ax is close to xo so that xo is a good initial guess to start the iterations. Once 
the exact new solution is found, the procedure can be repeated to trace the complete 
solution branch that plots the variation of x with respect to X starting from xq and Ao. This 
is the basic idea behind continuation methods. 

In case the function h(x, A.) is simple and there is a one to one correspondence 
between x and X, the above method can be applied with no problems. Complications 
however arise when there are multiple solutions x that satisfy h(x, A,) = 0 for a particular 
value of the parameter X. The case where all these solutions lie on the same solution 
branch indicates the presence of turning points with respect to X. This is illustrated in 
Figure 3.1 where there are two turning points with three steady states between A, = a and X 
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- c. These three steady states are shown for X = b in the figure. The two turning points are 


also labeled. 


Turning Points in a Solution Diagram 



Tracing the solution diagram around the turning points is complicated as when the 
parameter value is increased beyond the turning point, the previous converged solution is 
no longer a good starting guess for the Newton Raphson (or other) equation solver. Even 
if convergence is obtained, the middle portion of the solution diagram will be completely 
missed. This is illustrated in Figure 3.2. Two popular methods, namely, dynamic 
reparameterization and arc-length continuation are employed in the literature to trace the 
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solution branch in a reliable manner around a turning point. These are discussed in the 
following sub-sections. 


Branch Tracking Problems Near Turning Points 



Figure 3 .2: Illustration of problems in branch tracking near turning points 


3.1.1 Dynamic Reparameterization 

Dynamic reparameterization is in some sense, a “trick” method that utilizes the 
fact that a turning point with respect to X is generally not one with respect to at least one 
of the variables xj in x. Thus rather than choosing a value of X and then solving for x, we 
choose a value of xj and then solve for X and the remaining x variables. In other words, X 
is replaced with xj as the continuation parameter around turning points. The concept is 
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illustrated in Figure 3.3 for a scalar x variable. Note that at a turning point dx/dA, = co. 
Choosing the variable xj with the smallest rate of change dxj/dA, as the new continuation 
variable is generally a good way for solution branch tracking tracking. 


Dynamic Reparameterization 



Figure 3 .3 : Illustration of dynamic reparameterization for branch 

tracking around turning points 

Mathematically, let Xj.,A+ denote x with its element replaced by X ic 
Xj-.A+ = [Xi X2. . .Xj.i X Xj+i . . . •Xn]'^ 

Then if Xj is used as the continuation parameter around a turning point, a value of xj close 
to the previous point is chosen and the new value for Xj._x+ is found via Newton Raphson 
iterations. The Newton Raphson update formula is 

Xj.,x+''^' = Xj.,x+‘' + t.[Jj.x+]'* .h(xj.,x+) (3 -Sa) 
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where t is the step-correction between 0 and 1 and Jj.,x+ is the Jacobian of h with respect 
to Xj.,x+ and is obtained as 

'dh,ldx^ ••• dh,ldxj_, dh,ldX d\ldxj^, 


’J-M 


5/z, 1 8x^ 


dhj/dx, ••• dhj/dxj_, dhj/dX dhj/dxj^, ••• dhjldx„ 


[dhjdx, - dhjdxj_, dhjdX dhjdXj^, ■■■ dhjdx„ 


Once the turning point is crossed, the standard method of choosing a value of X 
and then solving for x is used. The Newton Raphson iteration formula in this case is 

x'^'■' = x'^^^t.[J]■^h(x^X'^) (3.2b) 

Here the Jacobian J is defined as 


J = 


d\/dx^ ■■■ d\ I dx, d\ I dx„ 


dhj / 5x, • ■ • dhj / dXj ■■■ dhj / 9x„ 


dhjdx^ ••• dhJdXj ••• dh„ldx„ 


Note that the Jacobian Jj. x+ is the same as J except that the column in J is replaced by 
dh/d?w. Limits on the slope dx/d X may be used to trigger reparametrization as the solution 
branch is tracked. Thus if dx/d X > tolerance 1, xj is used as the continuation and once the 
turning point is crossed and dx/d X < tolerance2, we continue with respect to X as the 
continuation parameter. Note that sign of the slope dx/dX determines if a positive or 
negative incremental change is made to X. A similar argument also holds for other 
continuation variables. 

It is quite clear from the discussion that the idea of dynamic reparameterization is 
very simple and easily programmable. We have therefore used it in this work on steady 
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state reactive distillation modeling. For the sake of completeness, a brief overview of the 
other popular method, arc-length continuation, is provided in the next sub-section. 


3.1.2 Arc Length Continuation 

In this case, the arc-length j is used to parameterize the solution branch. For the 
homotopy function h(x, X), no change in the function value occurs if we move exactly on 
the solution branch. In other words 


dh = hx.dx + hx.dA, = 0 (3.3) 

along the solution branch with hx = dh/dx and hx = dh/dX. Dividing by the differential 
arc-length ds and putting the equation in matrix form, we get 


[h. K] 


= 0 


(3.4a) 


where 


X = dx/ds 
and 

X = dX/ds 

For an n-element x vector, there are n+1 total variables (xandl) in equation 3.4a and 
only n relations. An additional relation is defined by the definition of arc-length from 
P 5 l;hagoras theorem as 

x2+X'=l (3.4b) 

Equations 3.4a and b form a square set of n+1 equations for n+1 variables. However, the 
arc-length definition (3.4b) is non-linear so that iterative methods may be necessary for 
solving the equations. The difficulty is however mitigated if we choose an arbitrary value 


25 



fori , solve for [x ly and then normalize [x 1]^ to unit length to satisfy the arc-length 
definition constraint. For a chosen value ofl = 1, from equation 3.4a, we get 


K 


X 


'O' 

0 

1 _ 



_ 1 _ 


(3.5a) 


and obtain 


X 


\ 

1 

-1 

'O' 

1 _ 

11 

0 

1 _ 


_ 1 _ 


(3.5b) 


The vector [ x 1 ]u^ gives the direction of the rate of change of x and A with respect to the 
arc length as we proceed along the solution branch. The subscript u is used to emphasize 
that the vector is unsealed. For standardization, this direction vector is normalized to a 
unit vector. The tangent direction unit vector at a point on the solution branch is then 
obtained as 

[X If =[x l]//||[x (3.5c) 


This tangential direction is also referred to as the Euler Predictor. If we proceed along 
this Euler Predictor direction with small enough arc-lengths, the complete solution branch 
can be tracked. In other words, if a small step As is taken along the Euler predictor 
direction, a new point on the solution branch is obtained as 




As 


In order to avoid accumulation of errors as the branch is tracked, the Newton Raphson 

method is used with [x AjnexJ as the initial guess for the iteration variables to ensure 
accurate tracking of the branch. In the above formulation, turning points are still a 
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problem as near a turning point, x obtained from equation 3.5 becomes very large as 
dX/ds = 0. As in dynamic reparameterization, a different slope Xj should be set to 1 

instead of 1 in equation 3.5(b) to calculate the tangential direction. Such 
reparametrization would ensure that the branch tracking progresses smoothly. 

This completes a brief overview of popular homotopy continuation algorithms in 
the literature. The interested reader is referred to Seader (1987) for an excellent review of 
the method. The next section formulates function h(x, A.) for reactive distillation problem 
as used in this work. 

3.2 Reactive Distillation Problem Formulation 

Reactive distillation involves simultaneous reaction and separation in a single 
column. Figure 3.4 shows a schematic of a general column with top down tray 
numbering. In order to develop a steady state model for the column, it is assumed that the 
vapor and liquid streams leaving a tray are at equilibrium. The equilibrium tray model is 
thus used to model the trays and the partial reboiler. A partial condenser is assumed so 
that it can also be treated as an equilibrium tray. 

Figure 3.5 shows a tray with all the material and energy flows. In the 
nomenclature, the subscript i is used to denote the components while subscript j is used to 
denote the tray number. The capital letters F, F and T denote feed, vapor and liquid 
material streams respectively. Corresponding small letters are used to denote component 
flows. The letter Q is used to denote energy streams. 
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Figure 3.4 A pictorial view of reactive distillation column 
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FIGURE 3.5: A schematic representation of an equilibrium stagej 

Of the various possible steady state model formulations, the Naphtali-Sandholm 
formulation has been used in this work. The reader is referred to Singh, (2004) for an 
excellent treatise on the Naphtali-Sandholm method for reactive distillation systems. In 
brief, at steady state the mole balances, equilibrium relations and heat balance must be 
satisfied for all the trays, reboiler and partial condenser. For and N-tray, C-component 
system with NRX number of reactions, if all the feed material streams, heat streams, 
stage pressures and side stream splits are specified, a square system of equations results 
with (N+2).(2C+1) equations and an equal number of variables. The governing MEH 
equations are 
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Mole Balances 


NRX 




/1=i 


Equilibrium Relations 


... (3.6) 


r} K I y 

I j 1,1 t,j 


kj 


E. 


Heat Balance 






•- V. , + 

'<J 


k = l 


= 0 


y V 

k,j-n 


(3.7) 


«, = K, 0 + ^;)I ‘u + \(1 + ■S,)I Vu - K E 

/=1 /=:I i = \ 

- K„ E ''u.i - E Sj = 0 


(3.8) 


In the above relations / = 1 to C and / = 1 to N+2 so that there are a total of 
(N+2).(2C+1) MEH equations. The above formulation assumes that the condenser and 
reboiler duties are specified since these are tray heat specifications for tray 1 (condenser) 
and tray N+I (reboiler) respectively. In other ■words, the column degree of freedom is 2. 
Specifying the condenser duty and the reboiler duty is not a good idea as the two are 
related trough an overall energy balance around the column. Thus if values that cannot 
satisfy the energy balance are input, the algorithm would refuse to converge. Alternative 
specifications that make more sense and converge easily are certainly possible. This is 
done by replacing the enthalpy balances for the condenser and reboiler respectively by 
the appropriate specification equations. Table 3.1 lists the replacement equations for 
some common specifications. The replacement equations for a total condenser are also 
noted in Table 3.2. 
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Table 3.1 ; The replacement functions for H; and Hji/ 

Specification 

Replacement for Hj 

Replacement for H 

Distillate or Bottom rate 

c 

yv,,-n = o 

c 

11 

cc 

1 




X/,, -iLID)l^v,, =0 y v ,^ =0 






Condenser or reboiler temperature 


Component mole fraction in c c c c 

product y V,. , - (£ V,. , =0 S w ’ (Z ) 


Component flow rate m distillate 
or bottom 




Table 3.2: Replacement equations for total condenser 


Equation 

Partial 

Total 

Equilibrium 

c 

C 

equations 

II 

i ly 

y 

+ 



tj C lU 

^k,j 

p=\ 


k=\ 

C 

= 0 
c 

for first component 


*=1 

II 


for i = 1 : C 

/,! /,1 C 

I. 

ifc = l 

for i = 2: C 












In the Naphtali Sandholm method, the component vapor and liquid flows leaving 
a tray and the tray temperatures are used as iteration variables. The governing MEH 
equations are linearized and the Newton Raphson method is used to converge to a bona- 
fide steady state. The tridiagonal algorithm is used to arrive at reasonable initial guesses 
for the iteration variables. 


3.2.1 Homotopy Continuation With Respect to Catalyst Weight 

The mole balances in the governing equations contain a generation term due to 
reaction. In the limiting case when this term is set to 0, the problem reduces to ordinary 
distillation with no reaction. A solution to the problem is easily found using the Naphtali- 
Sandholm method. This is akin to setting the catalyst weight (or hold-up for homogenous 
reactions) to 0 so that no reaction occurs. The catalyst weight can now be slowly 
increased to reach the steady state solution at the design catalyst weight. In order to 
implement the continuation, a catalyst weight factor L is used in the M equations such 
that at L= 0, the reaction term vanishes and at L = 1 , the reaction term corresponds to the 
design catalyst weight. More formally, the M equation is modified as 


^ - f>,j 

NRX 

-L o r . =0 

Z- ^ j ,n 


.. (3.9) 


Continuation with respect L generates the solution diagram of interest that clearly shows 
the existence of steady state multiplicities. The Naphtali-Sandholm algorithm for a given 
value of L is used to obtain the converged solution. Around turning points, dynamic 
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reparameterization is used and the appropriate column of the Jacobian matrix is replaced 
by dh/dL. 


dh nrx 

Note that — = - V u, r, 

dL it 


■j." 


(3.10) 


3.2.2 Homotopy Continuation About the Design Steady State 

Reactive distillation columns are designed for high purity products with high 
reaction conversion. Given that multiple steady states are routine in RD systems, 
appropriate control strategies must be implemented to ensure that the product purity and 
reaction conversion is maintained for all throughputs and expected feed composition 
changes. This requires a systematic study of the different input-output relationships in the 
column. 

For this purpose, continuation about the design steady state with respect to one of 
the specification variables such as reflux ratio, reboiler duty, tray temperature or tray 
composition is used to study the input-output relationships around the design steady state. 
In real-time operation appropriate control loops are used to maintain the column at these 
specifications. The specifications that result in simple input-output relationships with no 
multiplicities while maintaining product purity and reaction conversion are the ones that 
should be implemented. The corresponding control structure would be most robust in 
regulating the column for different production rates and disturbances such as feed 
composition change. Homotopy continuation techniques are necessary as the input-output 
relationships can be highly non-linear with output / input multiplicities. Output 
multiplicities occur when there are more than one values for the output variables such as 
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reaction conversion, product purity and tray temperature / composition profiles for a 
particular specification such as the reboiler duty, reflux rate and reflux ratio. Input 
multiplicities occur when there are multiple input specifications that give the same value 
for an output variable. 

It is noted that continuation about a design steady state is easier to implement than 
the catalyst weight homotopy continuation described earlier. This is because a standard 
RD simulator can be used to obtain the solution diagram as one of the specification 
variables is changed. It is also highlighted that the various choices of specifications in the 
Naphtali-Sandholm method would ensure accurate branch tracking around turning points 
using dynamic reparameterization. 

This completes a brief overview of the homotopy continuation method and its 
formulation for reactive distillation as applied in this work. The next two chapters 
describe, respectively, the application of catalyst weight continuation to methyl acetate 
and MTBE RD systems and control structure screening for an MTBE column using 
continuation with respect to a column specification about the design steady state. 


34 



Chapter 4 

Homotopy Continuation with respect to Catalyst Weight 

In this chapter, homotopy continuation with respect to catalyst weight is applied 
to two RD systems of commercial interest, namely, MTBE and methyl acetate. In the 
MTBE system, for a specified reflux ratio and bottoms flow rate, upto five steady states 
are detected in the solution diagram on two distinct branches. Two of these steady states 
disappear in the limit of reaction equilibrium on reactive stages at high catalyst weight. 
For the methyl acetate system, at fixed reflux rate and reboiler duty, three distinct steady 
states are seen in the solution diagram on two separate branches. In case the column 
specification is changed to fixed reflux ratio and reboiler duty, the solution diagram 
shows only a single steady state for all catalyst weights. Since each specification set 
corresponds to a particular column control strategy, these results show that the number of 
■ steady states can be affected by the control structure used and the catalyst weight per tray. 
The same must be chosen appropriately to minimize the steady state multiplicities for 
robust column operation and control. 

4.1 Case Study I: Methyl Tertiary Butyl Ether (MTBE) 

MTBE is an important industrial chemical used as a key octane enhancing 
ingredient in unleaded gasoline. Isobutene reacts with methanol to produce the desired 
product, MTBE, in the presence of an inert component, n-butane. The reaction is as 
follows 
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(CH, \C = CH, + CH^OH ^=± (CH, \ COCH, 

Isobutene MeOH MTBE 

An activity- based rate model expression is used to describe the kinetics of MTBE 
synthesis catalyzed by an ion exchange resin. The kinetic expression from Seader and 
Henley (1998) is used after replacing the component mole fractions with component 
activities. The modified rate expression for the forward reaction is 

^forward ~~ 3.67x10'^ exp(-92440/i?r)a/5 

The corresponding backward rate law is 

^backward “ 2.67 X lO'^ exp(-134,454/ RT)a^jgg I cd ueon 

where, 

= Yih 

The units of r are in moles per second per equivalent of acid groups. The catalyst 
is a strong ion exchange resin with 4.9 equivalents of acid groups per kilogram of 
catalyst. T is in kelvins and is the liquid mole fraction. The liquid-phase activity 

coefficients yi are modeled using the Wilson equation. The Wilson binary interaction 
parameters are taken from Doherty (2002) and are tabulated in Table 4.1. 


Table 4.1 ; Wilson interaction parameters for MTBE system 


Wilson Parameters 

IB 

MeOH 

MTBE 

NB 


IB 

0.0 

-0.7420 

0-2413 

0.0729 

MeOH 

0.7420 

0.0 

0.9833 

0.8149 

MTBE 

-’0.2413 

-0.9833 

0.0 

-0.1684 

NB 

-0.0729 

-0.8149 

0.1684 

0.0 

Kj 

IB 

0.0 

-85.5447 

30.2477 

0.0 

MeOH 

-1296.719 

0.0 

-746.3971 

-1149.280 

MTBE 

-136.6574 

204.5029 

0.0 

0.0 

NB 

0.0 1 

-192.4019 

0.0 

0.0 
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4.1.1 Column Configuration: 

The reactive distillation column shown in Figure 4.1 is used as the basis for this 
study. A similar column has been studied by Wang et al. (2003). The column consists of 
1 5 trays, a total condenser, and a partial reboiler. There are three zones in the column, a 
rectification zone (tray 1 and 2), a stripping zone (trays 11-15) and a reactive zone (trays 
3- 10). Mixed-butenes vapor at 350 K comprising 36% isobutylene and 64% inert n- 
butane by mole is fed to the column on tray 10. Pure methanol liquid at 320 K is fed on 
tray 9. The flow rates of methanol and mixed-butenes are nearly stoichiometric at 712.8 
and 1969.2 kmol/hr, respectively. The reflux ratio and bottoms flow rate are specified as 
7 and 640.8 kmol/hr, respectively. Pressure throughout the column is assumed to be 11 
atmospheres. High purity MTBE leaves in the bottoms stream while the unreacted 
isobutene and methanol and the n-butane inert are removed in the distillate stream. 




Figure 4. 1 : Reactive distillation column configuration for MTBE synthesis 
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4.1.2 Multiple Steady States by Continuation with Respect to Catalyst Weight: 

Homotopy Continuation with respect to catalyst weight is used to obtain the 
complete solution diagram. The reflux ratio and bottoms flow rate are assumed fixed 
as in Figure 4.1. In order to easily see the steady state multiplicities, the conversion 
and bottoms MTBE purity are plotted against L, the catalyst weight factor in Figure 
4.2. Two distinct branches are seen in the solution diagram for both purity and the 
conversion. Branch 1 is obtained by starting from the non-reactive solution and then 
introducing reaction by slowly increasing L. A starting point on branch 2 is obtained 
by solving the governing steady state RD equations using the traditional Naphtali- 
Sandholm method for the design catalyst weight (L=l). Continuation with respect to 
L about this steady state solution traces the second branch. 

Branch 1 has a single turning point at L = 0.705 while Branch 2 has two 
turning points at L = 1 .068 and L = 0.293 respectively. These turning points divide 
the complete solution diagram into four distinct regions with different number of 
steady states. As the catalyst weight factor is increased, the number of solutions 
increase or decrease by 2 when a turning point is crossed. In the first region 
corresponding to low' catalyst weights (L < 0.293) only one steady state solution 
exists on Branch 2. On crossing the turning point corresponding to L = 0.293 on 
Branch 2, the number of steady states increases to three are seen in the second region 
for which 0.293 < L < 0.705. All these solutions lie on Branch 2. As the catalyst 
weight is further increased, two additional steady states corresponding to Branch 1 
are encountered so that there a total of five steady states in the third region for which 
0.705 < L < 1.068. Moving beyond the L = 1.068 turning point on Branch 2, two of 
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the steady state solutions on Branch 2 disappear so that only three steady state 
solutions remain in the fourth region with L > 1.068. Two of these solutions lie on 
Branch 1 while the remaining one lies on Branch 2. The number of steady state 
solutions thus increases from one at low catalyst weights to three at intermediate 
catalyst weights. Five steady states are seen for a narrow range of catalyst weights in 
the kinetically controlled regime which reduces back to three at high catalyst weights 
when reaction equilibrium is approached on the reactive stages. 

The design catalyst weight (L=l ), lies in region three so that there axe a total 
of five steady states. The five steady states are labeled in Figure 4.2. Steady states I 
and II, lie on solution Branch 1 while the remaining steady states. III, IV and V, lie on 
Branch 2. The distillate and bottoms composition and reaction conversion for these 
five steady states are tabulated in Table 4.2. Of all these steady states, the high 
conversion steady state (Steady State V) with a conversion of 90.28% and product 
purity of 99.92%, is the most important as process economics dictates that conversion 
be as high as possible to enhance product purity and minimize the recycle of 
unreacted reactants. 
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Table 4.2; The distillate and the bottoms composition and reaction conversion for five 
steady states MTBE system 


SS 

Distillate Composition 

Bottoms Composition 

Conversion 

IB 

MeOH 

MTBE 

NB 

IB 

MeOH 

MTBE 

NB 

I 

0.2791 

0.0110 

0.0015 

0.7084 

0.0000 

0.8243 

0.1757 

0.0000 

15.95% 

II 

0.0970 

0.0156 

0.0006 

0.8868 

0.0000 

0.1515 

0.8485 

0,0000 

76.67 % 

III 

0.1008 

0.1138 

0.0001 

0.7853 

0.0062 

0.0257 

0.8193 

0.1488 

74.02 % 

IV 

0.0921 

0.1141 

0.0001 

0.7937 

0.0054 

0.0001 ^ 

0.8471 

0.1474 

76.53 % 

V 

0.0428 

0.172 

0.0003 

0.9397 

0.0000 

0.0001 

0.9992 ^ 

0.0007 

90.28 % 


In order to better understand these steady states, Figure 4.3 plots the product 
formation rate on the reactive trays for the five steady states. The temperature and 
components composition profiles are also shown in the Figure 4.4 and Figure 4.5 
respectively. Since the kinetic expression is quite complex with the methanol activity 
being in the denominator, it is hard to interpret the reasons for the very different 
product formation rates for the five steady states. 

For the low conversion steady state (SS I), a substantial temperature increase 
of more than 30 K occurs on the reactive trays (4-11). Since the activation energy of 
the backward reaction is larger than that for the forward reaction, the backward rate 
constant increases about 50 times while the forward rate constant only increases by 
about 15 times. Also, since for the backward reaction, the reaction order with respect 
to methanol is -2, the aa/aa^ term remains more than 1 as we move to tray 1 1 . Thus 
even though substantial forward reaction occurs on reactive stages 4 to 7, the net 
production is nearly zero as the backward reaction is favored on trays 8 to 10 due to 
larger backward reaction rate. The increase in temperature for SS 2 is less severe so 
that the backward reaction on trays 8-1 1 is relatively lesser than the forward reaction 
on trays 4-8. An intermediate conversion of 76.7% is therefore obtained. 
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The next two steady states, SS III and IV are very similar to each other in that 
a negligible change in temperature occurs on the reactive trays 4-11. The product 
formation profile is therefore flat. A net formation of MTBE occurs on all the reactive 
stages. Since the tray temperatures are lower compared to all other steady states, the 
product formation rates are also lower giving again an intermediate conversion of 
74% and 76.5%, respectively. For the high conversion steady state, SS V, higher 
product generation occurs for stages 4 to 10 with net consumption occurring on stage 
1 1 . The consumption on tray 1 1 is due to the substantial temperature rise across the 
reactive stages. Also, the low methanol mole fraction on tray 1 1 causes to be 
nearly 4 leading to a net backward reaction on tray 1 1 . The overall effect of higher 
product generation on all the reactive trays except tray 1 1 is a net conversion of about 
90%. The Activity ratios, forward and backward reaction rate are tabulated in Table 
4.3. 

For the MTBE system, several researchers have reported multiple steady 
states. Jacobs and Kfrishna (1993) found a high and a low conversion steady state. 
Chen et al. (2002) use catalyst weight homotopy continuation to report another steady 
state with intermediate conversion. Of the total three steady states, two lie on one 
solution branch and the last one on another solution branch. Sneesby et al (1998) 
trace the solution diagram by varying the reflux ratio at constant reboiler duty. They 
report a narrow zone of operating conditions with five steady states. Oiir results 
therefore are in qualitative agreement with these different literature reports. 
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Table 4.3. Activity ratio, forward reaction rate, back ward reaction rate and total reaction rates for five steady states of 
MTBE system 
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The results show that the interaction of reaction and separation can lead to 
enough non-linearity in the system so that steady state multiplicities occur. Note that 
when the reaction conversion is small, only a single steady state exists. In contrast, at 
higher conversions, three or five steady states exist. This shows that the reaction term 
in the governing equations is a major source of non-linearity. In practice, it is 
desirable to operate the column at high conversions such that the number of steady 
states in an open loop column is as low as possible. The task of regulating a column 
using appropriate control loops to maintain the product purity and conversion in the 
face of disturbances is easiest in case there are no steady state multiplicities and 
becomes progressively difficult as the open loop multiplicities increase. In the context 
of the MTBE column studied, it would therefore be desirable to provide a large 
amount of catalyst since there are only three steady states in contrast to five steady 
states at intermediate catalyst weights. Designing a column for still lower catalyst 
weights is not acceptable as the conversion and product purity become unacceptable. 
This is another example where operability and control considerations necessitate 
over-design. 



4.2 Case Study II: Methyl Acetate 

Methyl acetate is produced by the esterification of acetic acid with methanol 
in the presence of sulphuric acid or an acidic ion exchange resin. 


CH.COOH + CH^OH ^==± CH^COOCH^ + H^O 
HOAc MeOH MeOAc 

An activity based reaction rate model for heterogeneous catalysis, developed by Song 
et al ( 1 998), is used in this work. The rate expression is 


^catK 


R 


^HOAc ^MeOH 


^MeOAc 


■ K. 


MeOAc 


eq 


(l + + KMeOnaMeOH ^MeOAc^MeOAc ) 


where a denotes the component activities and Meat is the catalyst weight on a tray. 
The reaction rate constant, ki and equilibrium constant, Keq, are modeled as 

^- 6287 . 7 '^ 


k, = 6.942 X 10'° exp 
= 2.32 exp 


V T 
^ 782 . 98 ^ 


V 




Constant values for the adsorption equilibrium constants are used with 

=3.18, = 0.82,7f^,^,, =0.82, =10.5 


a. = y,x, 


In the above expressions, the reaction rate constant is in mole/ g cat / h, the catalyst 
weight is in grams and the temperature is in K. The net reaction on a tray is thus in 
kmol/hr. 
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The unwanted dehydration of methanol to dimethyl ether (DME) and water is 
ignoied in this study, as is the case with most literature reports. Accurate prediction of 
the VLE is very important for simulation and design of the column. The VLE is 
affected by the dimerization of acetic acid in the vapor phase. This is modeled using 
Marek’s method (Marek (1955)). Since the liquid phase is highly non-ideal, the 
Wilson equation is used to model the liquid phase activity coefficients. The Wilson 
parameters for the components are tabulated in Table 4.4. 


Table 4.4: Wilson interaction parameters for methyl acetate system 


Wilson Parameters 

HOAc 

MeOH 

Water 

MeOAc 


HOAc 

0.0 

-0.2583 

-1.1582 

0.3275 

MeOH 

0.2583 

0.0 

-0.8999 

0.5859 

Water 

1.1582 

0.8999 

0.0 

1.4858 

MeOAc 

-0.3275 

-0.5859 

-1.4858 

0.0 

K, 

HOAc ' 

0.0 

-1275.9 

-119.5 

-565.2 

MeOH 

275.6 

0.0 

-54 

-409.3 

Water 

-331.2 

-236.3 

0.0 

-965.4 

MeOAc 

350.4 

-15.7 

-325.0 

0.0 


4.1,1 Column Configuration: 

The reactive distillation column shown in Figure 4.6 is used as the basis for 
this study. A similar column has been studied by Al-Arfaz and Luyben (2002). The 
column consists of 35 trays, a total condenser, and a partial reboiler. There are three 
zones in the column, a rectification zone (trays 1-7), a stripping zone (trays 26-35) 
and a reactive zone (trays 8-25). Methyl acetate and water is formed between trays 8 
and 25. Pure acetic acid is fed to the column on stage 8 and pure methanol is fed on 
tray 25. Both feeds are saturated liquids and in exact stoichiometric proportion. The 
fresh feed rate is 300 kmol/hr. The pressure throughout the column is taken as 1 .25 

iPo A .1.4.8S73 - 
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Figure 4.6; Reactive distillation column configuration for Methyl acetate synthesis 
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atmospheies. At the design catalyst weight of 1000 kg / reactive tray, high purity 
methyl acetate leaves the column in the distillate while nearly pure water leaves in the 
bottoms. 


4.2.2 Multiple Steady States by Continuation with respect to Catalyst Weight: 

Homotopy continuation with respect to catalyst weight is used to study the 
steady state multiplicities in this RD system. A reflux rate specification of 581.25 
kmol/hr and a reboiler duty specification of 4.549 MW is used. Figure 4.7 plots the 
conversion versus the catalyst weight homotopy factor L. Two distinct branches are 
seen. These are labeled in the figure. Branch 1 is obtained by starting from the 
ordinary distillation solution at L = 0 and then introducing reaction by slowly 
increasing L. Branch 2 is traced by continuation about the reactive solution (L = 1) 
obtained using the conventional Naphtali-Sandholm simulator. In Branch 1, there are 
no turning points while Branch 2 has one turning point at L =0.15. The turning point 
divides the solution diagram into two regions, I and II, with one and three steady 
states respectively. For the design catalyst weight at L=l, there are a total of three 
steady states viz a high conversion, an intermediate conversion and a low conversion 
steady state. The distillate and bottoms compositions and the reaction conversion for 
these steady states are tabulated in Table 4.5. As in the MTBE case, it is the high 
conversion steady state that is of practical importance. 




Figure 4.7; N4ultiple steady states in methyl acetate system for fixed reflux rate and 

reboiler duty 


Table 4.5; The distillate and bottoms composition and reaction conversion for three 
steady states 


ss 

Distillate Composition 

Bottoms Composition 

Conversion 

HOAc 

MeOH 

Water 

MeOAc 

HOAc 

MeOH 

Water 

MeOAc 

High 

Conversion. 

0.0000 

0.0080 

0.0271 

0.9648 

0.0120 

0.0036 

0.9844 

0.0000 

98.83 % 

Intermediate 

Conversion 

0.0000 

0.0053 

0.0176 

0.9771 

0.1319 

0.1278 

0.7402 

0.0000 

85 . 10 % 

Low 

Conversion 

0.0000 

0.0041 

0.0001 

. 0.9829 

0.2271 

0.2249 

0.5427 

0.0054 

’ 70.08 % 
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In order to better understand the steady state multiplicities, Figure 4.8 plots 
the net rate of product formation on the trays for the three steady states. 


Reaction Rate 



Figure 4.8: Reaction rates for methyl acetate system 


The corresponding temperature profiles are also plotted in the Figure 4.9. 
Unlike the MTBE system, the temperature rise across the reactive section is only 10 
K, which is not substantial. 
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Figure 4.9: Temperature profiles and reaction rates for methyl acetate system 

The differences in the overall conversion are therefore largely determined by 
the composition profiles across the reactive section are shown in Figure 4.10. For the 
high conversion steady state, the mol fraction of acetic acid, a reactant, is 
substantially higher on trays 16 to 24 while the mol fraction of methyl acetate, a 
product, is substantially lower than the other two steady states. This leads to higher 
product formation rates resulting in nearly complete conversion. For the intermediate 
and low conversion steady states, the temperature profiles for are very similar so that 
the rate and equilibrium constants are the same. 
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Figure 4.10^ Cornposition profile for three stea.dy states of 
methyl acetate system 
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An examination of the net formation rate profile indicates that the net 
geneiation rates are nearly the same for all the reactive trays except trays 24 and 25. 
The leason is that for the intermediate conversion steady state, the unreacted 
methanol composition is significantly higher than for the low conversion steady state. 
This causes a larger forward reaction leading to nearly 85% conversion. 

As an aside, in order to study the effect of the column specification variables 
on steady state multiplicities, the specification set was changed to fixed reflux ratio of 
1.875 and reboiler duty of 4.549 MW. These specification values are equal to the 
corresponding high conversion steady state values obtained previously. The solution 
diagram corresponding to this specification set is shown in Figure 4.1 1. 


Fixed reflux ratio 



Figure 4.1 1; Multiple steady states in methyl acetate system for fixed reflux ratio and 

reboiler duty 




No steady state multiplicities are seen. It may be that another branch does 
actually exist but has not been traced because a starting point on that branch has not 
been obtained. We have used the Naphtali Sandholm RD simulator with several 
initial guesses for the iteration variables. It however always converges onto the 
branch shown in Figure 4.11. The results thus show that the column specification 
used can impact the number and also the regions of steady state multiplicities. Given 
that each column specification set corresponds to a particular column operating 
policy, for example, a fixed reflux rate or a fixed reflux ratio, the results show that the 
column operating policy determines the number of steady state multiplicities. This 
has important implications on the dynamic operation and control of an RD column. 
This aspect will be demonstrated in detail in the next Chapter using the MTBE 
column as an example. 

It is noted that in both the MTBE and methyl acetate case studies, in case multiple 
steady states occur, the reactive solution obtained using the NS method lies on a 
branch that is not connected to the branch that starts from the ordinary distillation 
solution. Thus homotopy continuation from a non-reactive solution provides 
additional steady state solutions that otherwise would be missed using conventional 
algorithms such as the NS method and the Inside-outside method. The synergy 
' between the computationally intensive homotopy continuation method and the 
traditional NS method in detecting multiple steady states in an RD colunrn is thus 
obvious. 
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Chapter 5 

Homotopy Continuation for Control Structure Synthesis: 

An MTBE Case Study 


An RD column must be operated so that the product purity and reaction 
conversion are maintained close to their design values for major disturbances entering the 
column. A robust control system needs to put in place that can regulate the column for 
anticipated production rate changes and changes in the feedstock composition, the 
principal disturbances into a column. In the design of such a control system, the selection 
of the control structure is the most crucial decision. 

A control structure refers to the number of control loops and the specific input- 
output pairing used in the loops. Potential input variables are the reflux rate, reflux ratio, 
reboiler duty, reboil ratio, distillate rate and bottoms rate. Potential output variables are 
easily measurable variables such as tray / stream temperatures and compositions. There 
are thus several possible input (manipulated) variables and output (controlled) variables 
even in a simple RD column. The various permutations and combinations lead to a large 
number of possible control structures from which a small set of “good” control structures 
must be chosen. Evidently, the key to successful colunui operation is this screening of 
control structures to zero in on the most appropriate controlled variables and effective 
handles to manipulate them so that by maintaining the controlled variables at their set- 
points, the purity and conversion of the column remain near the design specifications for 
the primary disturbances. Indeed, for a good control structure, other control system 
design decisions such as the choice between sophisticated versus simple control 
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algorithms becomes obvious. For a good control structure, simple PI control should 
provide the adequate column regulation. On the other hand, no amount of sophistication 
in the control algorithm can compensate for an inherently poor control structure. 

A good control structure is one that rejects disturbances effectively. In order to do 
so, controlled variables that are sensitive to the occurrence of the primary disturbances 
should be chosen so that timely control action can be initiated to regulate the column. 
Also, the manipulated variables used for control should affect the controlled variables in 
a substantial and easily predictable way. Such a choice assures enough “stick” for column 
regulation even for large disturbances. At steady state, linear or nearly linear input-output 
relationships are certainly desirable so that the process gain does not change appreciably 
and simple control algorithms suffice. In case multiple control loops are used, the input 
output pairings should not cause stability problems due to interaction between the loops. 
Last but not least, since sensor measurements are never exact, the column performance 
metrics such as reaction conversion and product purity should not be very sensitive to the 
set-point chosen for the controlled variables. In other words, the column performance 
should be robust to typical measurement errors or errors in the set-point inputs. 

Systematic steady state analyses can be conducted to reveal the control structure/s 
that provides effective column regulation for the primary disturbances in a column. These 
analyses include sensitivity studies for obtaining output variables that are sensitive to 
disturbances as well as manipulated variables. Sensitivity is an inherent property of a 
variable. Thus, a variable sensitive to disturbance will also be sensitive to manipulated 
variable. Note that each column specification set corresponds to a particular control 
structure. 
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Figure 5.1 shows the control structures corresponding to two different 
specification sets 



CS1 : Control structure for fixed reflux ratio and 
reboiler duty as specifications 



CS2: Control structure for fixed reflux ratio and tray 
temperature as specifications 


Figure 5.1: Control structures corresponding to two different specification sets 
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In reactive distillation systems, the coupling of reaction and separation in the 
same column leads to a highly non-linear system. The steady state input-output 
relationships aie complex and multiplicities are known to occur. We distinguish between 
two types of multiplicities; output multiplicity and input multiplicity. The former refers to 
multiple output values for the same input specification while the latter refers to multiple 
input specifications giving the same output. This was illustrated in Figure 1.4. From a 
control perspective, input variables are potential manipulated variables such as the 
reboiler duty or reflux ratio while output variables are potential controlled variables such 
as tray temperatures / compositions etc. The choice of the control structure substantially 
affects the number of steady state multiplicities and the size of the operating space in 
which these multiplicities occur. In other words, the complexity of the input-output 
relationships depends on the control structure. The control structure chosen must avoid 
output multiplicities to ensure that the column does not drift from say a high conversion 
steady state to a low conversion steady state even though there is no apparent change in 
the inputs to the column. Input multiplicity must be avoided as the sign of the gain 
between the input and output variable changes depending on the column operating 
condition. Such changes in the sign of the gain can easily lead to dynamic instability 
problems. Also, if the output variable is used in a control loop, the range of feasible set- 
point specifications becomes limited so that measurement biases can drive the control 
system to seek a infeasible steady state. 

Systematic steady state analyses can be conducted to formulate good control 
structures that would regulate the column in the presence of disturbances. Two distinct 
approaches are possible to control structure design. In the first, all reasonable candidate 
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control structures are first enumerated and then accepted or rejected depending on their 
ability to maintain the column close to the design purity and conversion for the primary 
disturbances. Structures with linear input-output relationships and no steady state 
multiplicities are preferred. Alternatively, in the second approach, the effect of varying 
the input handles about their base case values, on the column outputs is first studied to 
identify potential input-output pairings. Pairings that are sensitive avoid multiplicities and 
give nearly linear input output relationships are short listed. Next, control loops are 
systematically added and the column steady state response to the primary disturbances 
evaluated until we arrive at the simplest control structure that provides adequate 
regulation in the presence of disturbances. 

This chapter demonstrates the application of the latter approach to the synthesis of 
a “good” control structure for a double feed MTBE RD column. The preferred control 
structure operates at fixed reflux ratio and uses a temperature inferential control loop for 
maintaining the product purity by manipulating the reboiler duty. A reactive tray 
composition is controlled for balancing the fresh feed into the column as per the reaction 
stoichiometry by manipulating a fresh feed rate. The isobutene fresh feed or the methanol 
fresh feed can be used as the manipulated variable in the composition loop. The tray 
location for composition measurement must be carefully chosen to avoid stability 
problems due to loop interactions. 

5.1 The MTBE Reactive Distillation Column 

The double feed column described earlier in Chapter 4 is used for the control 
structure synthesis case study. A schematic of the column and the base case operating 
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conditions are shown in Figure 4.1. A similar column has been studied by Wang et al 

(2003) for control strategy of reactive distillation for MTBE synthesis. Reaction kinetics 
used by same as follows: 


r = MJi^ 


a 


IB 


L ^MeOH 


^MTBE 

K rP' 

MeOH 


A:;- — 3.67x lO^^e moll{sequiv) 


= exp(- 16.33 + 6820 /r) and <7. = y.x^ 


Where represents the reaction equilibrium constant and T is temperature in K. 

Corresponding to this reaction kinetics homotopy continuation with respect to catalyst 
weight is applied shown in Figure 5.2. 


Multiple steady states in MTBE system 



Figure 5.2: Multiple steady states in MTBE system for purity and conversion 
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This solution diagram also shows five steady states for a unique catalyst weight. 

1 his diagram is ditferent from earlier one which also shows five steady state solutions. 
The difference is attributes to the different kinetic parameters as described in earlier 
chapter 4. The catalyst weight must be chosen appropriately to minimize the steady state 
multiplicities for robust column operation and control. Therefore, 1800 kg catalyst (L = 

1 .2) is chosen at which only three steady states are present with high conversion purity. 

1 he column operates at a high reflux ratio of 7 with a base case conversion of 
90% with >99% pure MTBE leaving from the bottoms. A control system that maintains 
the product purity and reaction conversion near the design values is desired. It should 
provide smooth transitions for production rate changes and compensate for any changes 
in the butene feed composition. In order to synthesize such a control structure the effect 
of the various potential manipulated variables on the column tray temperatures and 
compositions is studied. 

5.2 Effect of Column Inputs on Output Variables 

The potential input variables or valves that can be manipulated to regulate the 
column are the two fresh feeds, the two product streams, the reflux rate and the reboiler 
duty. It is also possible to maintain the flow ratios of two streams such as the reflux ratio 
or the reboil ratio. Manipulating the reflux rate or ratio for column regulation is not a 
good idea since the hold up (or catalyst weight) on each reactive stage in an RD. column 
is large so that the column response can be extremely sluggish. In contrast, manipulating 
the reboiler heat duty causes an immediate change in the vapor rate throughout the 
column. These dynamic considerations indicate that fixed reflux rate or reflux ratio 
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operating policies are preferable and the reboiler duty and feed flows can be used for 
regulating the separation and the reaction in the column. 

Now in ordei to decide between the fixed reflux rate vs fixed reflux ratio policy, 
Figure 5.3 plots the steady state purity and conversion as the reboiler duty is varied with 
the feed conditions at their base case values for the two operating policies. Output 
multiplicity is evident for both the fixed reflux rate and the fixed reflux ratio policies. It is 
however noted that in the latter case, the base case point is away from the zone of output 
multiplicity while output multiplicity occurs in the former even for the base case. The 
column can thus drift from a high conversion steady state to a low conversion steady state 
even as the reboiler duty and the feed conditions remain the same for a fixed reflux rate. 
The same cannot happen in case the reflux ratio is kept fixed. Clearly, for this column the 
fixed reflux ratio is preferable over the fixed reflux rate policy. 
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Having decided on maintaining a fixed reflux ratio, the response of the tray 
temperatures and compositions to changes in one of the remaining input variables is 
studied, figures 5.4 and 5.5 plot the steady state tray temperature and composition 
1 esponses to changes in the reboiler duty about the base case respectively. 


Temperature profiles 



Figure 5.4: Responses of tray temperatures to changes in reboiler duty at 
fixed reflux ratio 
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Figure 5.4; Responses of tray temperatures to changes in reboiler duty at 
fixed reflux ratio 
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Figure 5.5; Responses of tray composition to changes in reboiler duty at 
fixed reflux ratio 

All the tray temperatures and compositions, especially for trays 12 to 15 in the 
stripping section respond to changes in the reboiler duty. A small zone with output 
multiplicities is observed. The best tray location for a temperature sensor is the as the one 
that is most sensitive to changes in the reboiler duty. A plot of the dT/dQ also referred to 
as the tray temperature sensitivity with respect to reboiler duty profile in Figure 5.6, 
shows that tray 12 is the most sensitive to changes in the reboiler duty. The temperature 
or composition of a responsive tray may be controlled by mampulating the reboiler duty 
in order to regulate the separation, as is usually done in ordinary distillation. Temperature 
measurements are preferred over composition measurements due to their accuracy and 
quick response time. 
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Figure 5.6: Temperature sensitivity with respect to reboiler duty 

Next, the effect of varying the fresh feed flows on the column tray temperatures / 
compositions is evaluated. The tray temperatures / compositions for changes in the fresh 
methanol shown in figure 5.7 and 5.8 respectively whereas, for isobutene feeds at fixed 
reflux ratio and reboiler duty are plotted in Figures 5.9 and 5.10 respectively. Substantial 
input multiplicity in the temperature of reactive trays is evident from the plots. The 
reactive tray temperature should therefore not be controlled using the fresh feed as the 
process gain changes sign depending on the operating condition. Also the tiay 
temperature set point can be infeasible due to sensor bias. For exampl® a temperature set 
point of 370.03 K for tray 1 1 does not have a feasible steady state solution as it is above 
the maximum temperature as seen in the Figure 5.9. The reactive tray compositions in 
contrast are well behaved. 
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Figure 5.7: Tray temperatures for changes in MeOH feed flow 
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Figure 5.9; Tray temperatures for changes in Butenes feed flow 
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Figure 5.10: Tray compositions for changes in Butenes feed flow 
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A sensitive tray is needed for effective feedback control. Figure 5.11 plots the 
sensitivity of tray composition with respect to the fresh feed flow rates. Of the reactive 
trays, tray 4 is the most sensitive to changes in the fresh methanol feed while trays 7 and 
8 are the most sensitive to changes in the butene feed. 



Figure 5.11: Sensitivity of the tray compositions with respect to fresh feed flow 
rates 
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5.3 Control Structure Synthesis 

Wc sldit with the simplest column operation strategy of no control and then keep 
adding complexity in the form of control loops until the simplest control structure that 
achieves the desired control objective. The objective here is to maintaiTi the product 
purity and reaction conversion for the principal disturbances, namely production rate 
variation and feedstock composition. Operating the column at fixed reflux ratio and a 
fixed reboiler duty is the starting point. The fresh feeds are' assumed under flow control. 
The corresponding control structure is shown in Figure 5.1 and labeled as CSl. Figure 
5.12 plots the column purity and conversion as the fresh feed rates vary about the base 
case values. These feed rates may be varied by an operator or may be off due to the 
typically high errors involved in flow measurements. Even for relatively small changes in 
the fresh feed rates, the product purity becomes unacceptable. This is because one of the 
reactants becomes excess and exits with the product stream, reducing its purity. Clearly 
the feeds need to be kept in stoichiometric balance in order to avoid this problem. 


MTBE Purity and Conversion 






If the feeds are assumed to be in stoichiometric balance using some feedback 
arrangement, the column performance in terms of product purity and reaction conversion 
is still unacceptable as the production rate is changed ie more fresh feed is put into the 
column. This is shown in Figure 5.13. Thus operating the column at fixed reflux ratio and 
reboiler duty with stoichiometric feeds is not satisfactory. Appropriate control loops 
therefore need to be added. 


Stoiohiometrio feed 



Figure 5.13: Purity and conversion vs production rate change at fixed reflux ratio 
and reboiler duty with stoichiometric feeds 
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The rcboiler duty and fresh feeds need to be adjusted in a manner that is as yet not 
known. The previous section on input-output relationships showed that the stripping 
section tray tempeiatures are quite sensitive to changes in the reboiler duty. Accordingly, 
as in ordinary distillation, a tray temperature control loop that manipulates the reboiler 
duty to keep a tray temperature in the stripping section fixed, is added. The corresponding 
control structure is shown in Figure 5.1 and labeled as CS2. Since tray 12 is the most 
sensitive to the reboiler duty, it is chosen as the location for placing the temperature 
sensor. Figure 5.14 plots the conversion and MTBE purity for a fixed reflux ratio and a 
fixed tray 1 2 temperature as both the fresh feed rates are varied in ratio about the base 
case. Given that the fresh feed rates are nearly stoichiometric, the temperature controller 
ensures high conversion and purity for all feed rates, (the production rate handle). 


At stoichiometry feeds and fixed 12“' tray temperature 



Figure 5.14: Purity and conversion with production rate change at fixed reflux ratio and 
tray 12 temperature with stoichiometric feeds 
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As before, near stoichiometric fresh feeds flows are not possible due to large 
erioi in the flow measurements. Also, even if the feed flows are equal, composition 
changes in the flesh feeds can upset the stoichiometric balance. This can lead to 
unacceptable product purity as the unreacted excess reactant must exit the column. The 
deterioration in steady state column performance is illustrated in Figure 5.15 for feed 
flow imbalances and feed composition changes. 



Feed composition changes 



Figure 5.15: Purity and feed composition changes at fixed reflux ration and 
reboiler duty 
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It IS clear that a feedback control loop arrangement is necessary in order to 
maintain the feeds in stoichiometric balance. The temperature or composition in a 
leactive tiay may be controlled using one of the feed streams for this purpose. The 
previous section showed that due to input multiplicities, reactive tray temperatures should 
not be used in a control loop. The tray compositions, on the other hand, are potential 
controlled variables for balancing the fresh feeds. A tray composition control loop using 
one of the ftesh feeds as the manipulation handle is therefore added. The corresponding 
control structures are shown in Figure 5.16 are labeled as CSS and CS4. Note that since 
there are now two control loops, the reactive tray location for composition for 
temperature measurement and the fresh feed handle must be chosen to avoid instability 
due to interaction between the two control loops. 

The overall decision on which fresh feed and reactive tray composition to use in 
the control loop is now a compromise between tray sensitivity, low interaction and quick 
dynamic response time. The Niederlinski Index (NI), a measure of loop interaction, is 
tabulated in Table 5.1 for composition measurements in different reactive trays. A 
negative NI indicates that the two control loops will be closed loop unstable while a NI 
between 0.4 to 2 is considered acceptable. For the methanol fresh feed, tray 4 shows the 
highest sensitivity and has an acceptable NI of 0.51. The NI for other sensitive reactive 
trays is near zero or negative. For the butenes fresh feed, tray 1 1 is sensitive and has an 
acceptable NI of 0.46. Note that trays 6-9, though sensitive, have a negative NI and 
cannot be used for control. Thus a composition controller that either manipulates the 
methanol feed to track the tray 4 temperature or manipulates the butenes feed to track 
tray 1 1 temperature should balance the feed flows to satisfy the reaction stoichiometry. 
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Reflux Ratio 
Set Point 


Reflux Ratio 
Set Point 




CSS: MeOH as a production rate handle 


CS4: Butenes as a production rate handle 


Figure 5.16; Control structures with two control loops 


Table 5. 1 : Niederlinski Index corresponding to butenes and MeOH feed flows 


Tray Number 4 5 6 7 8 9 10 11 

Butenes 2.44 18.03 -1.55 -0.675 -0.513 -0.684 1.35 0.461 




In order to assess the stead, state eolumn performance of these control structures, 
the product purity and reaction conversion are plotted for production rate changes and 
feed composttion changes. Note that for CS3 the production rate handle is the fresh 
methanol feed while for CS4 the production rate handle is the fresh butene feed. Both the 
structures are able to maintain the column at high product purity and reaction conversion. 
Figure 5.17 and Figure 5.18 plot the column inputs (manipulated variables) with respect 
to the production rate for both control structures respectively. The relationships are linear 
for both the control structures so that simple PI controllers should provide good column 
regulation. Column input versus butene feed composition is plotted in Figure 5.19. Mild 
non-lineaiity in the input output relationships is observed. Simple PI controllers should 
again provide adequate column regulation for feed composition changes. 


89 



Bottoms flow rate (k mo I/hr) — Purity and conversion 










Bottom flow rate (kmol/hr) 


0.99 ■--- 
0.98 -■ 
0.97 


Feed composition changes 


4- Purity 
O Conversion 


It! 

> 0.94 


case pase = (|).36 


D.32 0.33 0.34 0.35 0.36 0.37 0.38 0.39 

Isobutene feed composition 


Feed composition changes 






-5 

^ase ca?e = 0.^6 

1 1 I 



1 

L - J 



u 

I I I 

+ : : 

I 1 1 

L 


0.32 0.33 0.34 0.35 0.36 0.37 0.38 0.39 0.4 0.41 

Isobutene composition In teed 

Figure 5.19: Column Performance with Production Rate Changes 
Production rate Handle: MeOH 






I he MTBE example clearly demonstrates how systematic steady state analyses 
can be used loi synthesizing control structures that ensure good column regulation in the 
presence of disturbances. These control structures need to be validated via dynamic 
simulations. This shall be taken up in the future work. It is highlighted that in a recent 
aiticle, Wang et al (200j) show the dynamic column performance of a control structure 
similar to CS3. The control structure is able to handle ±20% changes in the production 
rate without any problems. 

At this point, it is highlighted that several alternative control structures that are 
conceptually similar to CS3 are possible. These are shown in Figure 5.20 and labeled. Of 
the control structures shown (a) would provide the tightest column regulation. In (a), the 
two fresh feeds are ratioed and the reboiler duty is ratioed with the production rate 
handle. A reactive tray composition controller provides the ratio set-point to the former in 
a cascade arrangement for maintaining the stoichiometric balance. Similarly, the stripping 
section temperature controller provides the ratio set-point to the latter for maintaining 
product purity. The ratio controllers provide feed-forward compensation for production 
rate changes so that the master controllers in the cascade loops only have to make fine 
adjustments by correcting the ratio set-points. Tight control is thus possible. 
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Thio completes a demonstration of control structure synthesis for an MTBE RD 
column using systematic steady state analyses. The MTBE case study shows that 
maintaining the stoichiometric balance of the fresh feeds is the key to successful RD 
column operation. The results also show the substantial effect that a control stmcture has 
on steady state multiplicities. The next chapter draws conclusions from this work and 
provides directions for future work. 
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Chapter 6 

Conclusions and Future Work 


6.1 Conclusions 

Conventional methods such as Naphtali— Sandholm and Inside — Outside methods, 
depending on the initial guess can converge to any one of the multiple steady states, in 
case multiple steady states do exist. Few people obtained merly two steady states by 
changing initial guesses using these conventional methods. Therefore, these methods are 
not geared towards detecting steady state multiplicities. Homotopy continuation, a 
powerful tool for identifying steady state multiplicities, is applied to reactive distillation 
systems to assess the impact of these multiplicities on column design, operability and 
control. In the MTBE system, a reasonably large zone with five multiple steady states is 
obtained with the high conversion steady state corresponding to the design steady state. 
As the catalyst weight is further increased, two of these steady states disappear. Thus 
there were five steady states in the kinetically controlled RD regime in contrast to only 
three steady states in the reaction equilibrium controlled RD regime. For the methyl 
acetate column, three steady states in both the kinetically controlled and reaction 
equilibrium regimes are obtained at fixed reflux rate and reboiler duty. If the column- 
specification is changed to fixed reflux ratio and reboiler duty, only a single steady state 
solution is obtained for all catalyst weights. Therefore, the number of steady states is 
affected by the column specification used and the catalyst weight per tray. 

Homotopy continuation is used for control structure synthesis for MTBE system. 
The input-output relationships are systematically analyzed to synthesize robust column 
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operating structures that ensure high product purity and reaction conversion in the face of 
majoi distui bailees. Column specifications also found very important to avoid or 
minimized multiple steady state in control structure synthesis. In order to decide between 
the fixed reflux rate vs fixed reflux ratio as a column specification, Output multiplicity 
was evident for both the fixed reflux rate and the fixed reflux ratio specification. It was 
however noted that in the latter case, the base case point was away from the zone of 
output multiplicity while output multiplicity occurs in the former even for the base case. 
Clearly, for this column the fixed reflux ratio was preferable over the fixed reflux rate. 
Sensitivity analysis played an important role for obtaining the sensitive tray location for 
the temperature / composition sensors used in the control structure. The various 
permutations and combinations lead to a large number of possible control structures from 
which a small set of “good” control stmetures are chosen. The preferred control stmeture 
operates at fixed reflux ratio and uses a temperature inferential control loop for 
maintaining the product purity by manipulating the reboiler duty. A reactive tray 
composition is controlled for balancing the fresh feed into the column as per the reaction 
stoichiometry by manipulating a fresh feed rate. The isobutene fresh feed or the methanol 
fresh feed can be used as the manipulated variable in the composition loop. This case 
study shows that maintaining the stoichiometric balance of the fresh feeds is the key to 
successful RD column operation. The MTBE example clearly demonstrates how 
systematic steady state analyses can be used for synthesizing control structures that 
ensure good column regulation in the presence of disturbances 
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6.2 Future Work 


I his developed homotopy continuation algorithm is working manually. In future, 
this homotopy continuation algorithm should be in such a way that can give whole 
homotopy curve automatically. A pre defined step size is needed to trace the homotopy 
curve. A reliable step size estimation method that predicts the step length can reduce the 
computation time tremendously. George (1980), den Heijer and Rheinboldt (1981), 
Deuflhard (1979), and Wacker (1978), have proposed different procedures to es timate 
step size. These proposed control structures need to be validated via dynamic simulations. 
This shall be taken up in the future work. 
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